The inverse cascade and nonlinear alpha-effect in simulations of 
isotropic helical hydromagnetic turbulence 
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ABSTRACT 

A numerical model of isotropic homogeneous turbulence with helical forcing is investigated. 
The resulting flow, which is essentially the prototype of the a 2 dynamo of mean-field dynamo 
theory, produces strong dynamo action with an additional large scale field on the scale of the 
box (at wavenumber k = 1; forcing is at k = 5). This large scale field is nearly force- free and 
exceeds the equipartition value. As the magnetic Reynolds number R m increases, the saturation 
field strength and the growth rate of the dynamo increase. However, the time it takes to built 
up the large scale field from equipartition to its final super-equipartition value increases with 
magnetic Reynolds number. The large scale field generation can be identified as being due to 
nonlocal interactions originating from the forcing scale, which is characteristic of the a-effect. 
Both a and turbulent magnetic diffusivity rjt are determined simultaneously using numerical 
experiments where the mean-field is modified artificially. Both quantities are quenched in a In- 
dependent fashion. The evolution of the energy of the mean field matches that predicted by 
an a 2 dynamo model with similar a and n t quenchings. For this model an analytic solution is 
given which matches the results of the simulations. The simulations are numerically robust in 
that the shape of the spectrum at large scales is unchanged when changing the resolution from 
30 3 to 120 3 meshpoints, or when increasing the magnetic Prandtl number (viscosity /magnetic 
diffusivity) from 1 to 100. Increasing the forcing wavenumber to 30 (i.e. increasing the scale 
separation) makes the inverse cascade effect more pronounced, although it remains otherwise 
qualitatively unchanged. 



Subject headings: MHD - turbulence 

1. Introduction 

The generation of large scale magnetic fields 
from small scale turbulence is important in many 
astrophysical bodies (planets, stars, accretion 
discs and galaxies). Over many decades the auj- 
dynamo concept has been invoked to explain large 
scale magnetic field generation (Moffatt 1978, 
Parker 1979, Krause & Radler 1980). Over re- 
cent years however numerical simulations have 
become available that produce large scale fields 
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with appreciable magnetic energy, sometimes even 
exceeding the turbulent kinetic energy (e.g. Glatz- 
maier & Roberts 1995, Brandenburg et al. 1995, 
Ziegler & Riidiger 2000). Whether or not large 
scale field generation to such amplitudes is related 
to the a-effect remains debatable (e.g. Cattaneo 
& Hughes 1996, Brandenburg & Donner 1997). 

The a-effect is a key ingredient to many as- 
trophysical dynamo models. The purpose of this 
paper is, therefore, to study a simple system 
that is prototypical of the a-effect: homogeneous 
isotropic turbulence that lacks mirror symmetry. 
[Astrophysical dynamos often work in conjunction 
with shear, i.e. the w-effect: this case is stud- 
ied in a second paper (Brandenburg, Bigazzi, & 
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Subramanian 2000)]. An isotropic helical turbu- 
lent flow is accomplished by adopting a body force 
corresponding to plane polarized waves in ran- 
dom directions (but constant polarization) with 
wavelengths short compared with the size of the 
box. Since the seminal papers by Frisch et al. 
(1975) and Pouquet, Frisch, & Leorat (1976) we 
know that there should be an inverse cascade of 
magnetic helicity, which has also been demon- 
strated using direct numerical simulations (e.g. 
Meneguzzi, Frisch, & Pouquet 1981, Balsara & 
Pouquet 1999). However, to our knowledge there 
has never been a detailed study of the spatial mag- 
netic field patterns obtained from the inverse cas- 
cade, nor has there been a quantitative identifica- 
tion of the classical a-effcct in mean-field dynamo 
theory. Furthermore, the Reynolds and Prandtl 
number dependences of this process have not been 
fully explored yet. In the present paper we study 
models with strongly helical forcing at different 
Reynolds numbers. We also investigate some mod- 
els where the magnetic Prandtl number (viscos- 
ity/magnetic diffusivity) is increased from 1 to 
100. This may be important in connection with 
the galactic magnetic field, and there are some se- 
rious concerns that the inverse cascade may not 
be efficient at large magnetic Prandtl numbers. 

2. The model 

We consider a compressible isothermal gas with 
constant sound speed c s , constant dynamical vis- 
cosity p, constant magnetic diffusivity i], and con- 
stant magnetic permeability p - The governing 
equations for density p, velocity u, and magnetic 
vector potential A, are given by 



D lnp 
Dt 



= -V • it, 



(1) 



Dm 



-^Vln/)+^+^(V 2 «+iVV-«)+/, 



dA 
~dt 



= u x B rjnoJ, 



(2) 
(3) 



where D/Dt = d/dt + u ■ V is the advective 
derivative, B = V x A is the magnetic field, 
J = V x B/p is the current density, and / is 
a random forcing function. 

We use periodic boundary conditions in all 
three directions for all variables. This implies that 



the mass in the box is conserved, i.e. (p) = po, 
where po is the value of the initially uniform den- 
sity, and angular brackets denote volume averages. 
We adopt a forcing function / of the form 

f(x, t) = Rc{Nf k(t) exp[ik{t) ■ x + i<f>(t)}}, (4) 

where k(t) — (k x ,k y ,k z ) is a time dependent 
wavevector, x = (x,y,z) is position, and <p(t) 
with |(/)| < 7r is a random phase. On dimensional 
grounds the normalization factor is chosen to be 
TV = /oc s (fcc s /(5t) 1 / 2 , where /o is a nondimensional 
factor, k = \k\, and St is the length of the timestep. 
We focus on the case where |fe| is around fcf = 5, 
and select at each timestep randomly one of the 
350 possible vectors in 4.5 < |fc| < 5.5. We force 
the system with eigenfunctions of the curl opera- 
tor, 

k x (fc x e) — i\k\(k x e) 



fk 



2fc 2 yT - {k-ef/k 2 



(5) 



where e is an arbitrary unit vector needed in order 
to generate a vector k x e that is perpendicular 
to k. Note that |/^.| 2 = 1 and, in particular, 
ik x ffc = \k\ffc, so the helicity density of this 
forcing function satisfies 



/•Vx/ = |fe|/ 2 >0 



(6) 



at each point in space. We note that since the 
forcing function is like a delta- function in fc-space, 
this means that all points of / are correlated at 
any instant in time, but are different at the next 
timestep. Thus, the forcing function is delta- 
correlated in time (but the velocity is not). 

We adopt nondimensional quantities by mea- 
suring u in units of c s , x in units of 1/fci, where 
k\ is the smallest wave number in the box, which 
has a size of L = 2ir, density in units of po, and B 
is measured in units of y^poPo c s- This is equiva- 
lent to putting 



ki = po = po = 1. 



(7) 



In the following we always quote the mean kine- 
matic viscosity v = p/po, which is close to the 
actual kinematic viscosity /x/ p because the Mach 
numbers considered in the present paper are less 
than one. 

We advance the equations in time using a third 
order Runge-Kutta scheme and sixth order explicit 
centered derivatives in space. In all cases pre- 
sented we chose /o = 0.1, which yields rms Mach 
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numbers around 0.1-0.3, and peak values less than 
one. 

Our initial condition is lnp = u = 0, and A 
is a smoothed gaussian random field that is delta- 
correlated in space, so the initial magnetic energy 
spectrum is Eyi{k) ~ k 4 with a decline at high 
wavenumbers. 

3. Results 

All the runs are summarized in Table 1. The 
definition of various entries to the table are given 
below, together with an outline of the general be- 
havior of the solutions. 

After about 30 time units the rms velocity, w rms , 
reaches an approximate equilibrium amplitude of 
up to 0.3. (Since c s = 1, this is also the Mach num- 
ber.) This velocity corresponds to a turnover time 
of r — £{/u rms sa 4 time units, where £{ — 2w/k{ 
is the forcing scale. We note that the value of r 
is approximately equal to the value of the corre- 
lation time obtained from the temporal correla- 
tion function of the velocity. The value of u rms is 
somewhat smaller for smaller Reynolds number. 
The flow has strong positive helicity, as measured 
by the relative helicity (u> • u) / (oj lms u lms ) , which 
can be as large as 70% (or even larger when the 
Reynolds number is smaller). Here, us — V x u is 
the vorticity. 

The growth rate of the magnetic field is deter- 
mined as A = dln(B 2 ) 1 1 / 2 /dt, where the subscript 
'lin' refers to early times when the field is still weak 
on all scales. A hyphen in the table indicates that 
the run has been restarted from another run, so 
no data are available for the linear growth phase. 
Also given is the growth rate normalized with the 
turnover, r = lf/{u 2 )\{^. 

In order to assess the Reynolds number depen- 
dence of our results we have performed three runs 
with v = 7] = 0.002 (Run 1), v = r) = 0.005 
(Run 2), and v = t] = 0.01 (Run 3); see Table 1 
for a summary. In order to assess the dependence 
on magnetic Prandtl number we have additional 
runs with v/r] = 20 (Run 4) and 100 (Run 5), as 
well as one with v/r] = 0.1 (Run 7). These runs 
will be explained in detail in §3.5. In Runs 1-5 
and 7 the forcing wavenumber was around 5, but 
in Run 6 we increased it to 30 in order to study 
the properties of larger scale separation; see §3.5. 
The root-mean-square values of various quantities 



are reasonably well converged, as can be gauged 
by comparing Run 2 (60 3 meshpoints) with Run 21 
(30 3 meshpoints), which has the same values of rj 
and v. We return to a detailed discussion on the 
Reynolds number dependence in §3.6. 

In the table we give various magnetic Reynolds 
numbers: R m is based on the box size (= 2ir) and 
the velocity (u 2 )),^ at the time when the mag- 
netic field is saturated, i? mj n n is the same but dur- 
ing the linear growth phase (using (u 2 )\(^), R mt \ 
is based on the Taylor microscale (u 2 )\(^ / (u} 2 )\(^ 
and i? mj f rc is based on the forcing scale £{ and 
(it 2 )] 1 / 2 . The critical values of i? mj f or c for the onset 
of dynamo action are also given and are typically 
between 7 and 9. In all cases the onset for dy- 
namo action occurs for v/r] < 1, i.e. for magnetic 
Prandtl numbers less than unity. 

The evolution of magnetic and kinetic ener- 
gies (Em and Ek) is shown in Fig. 1. Note that 
Ek decreases after Em has reached its saturation 
value. (We note that even after saturation the 
field continues to grow somewhat, but this will 
be discussed in full detail in §3.6.) The rela- 
tive kinetic helicity changes only slightly before 
and after saturation. Both the growth rate and 
the saturation level of the magnetic field increase 
with increasing Reynolds number, and are likely to 
reach some asymptotic value at sufficiently large 
Reynolds number. 

The level of turbulence may be characterized by 
the ratio of the turbulent to the microscopic dif- 
fusion coefficient for a passive scalar, D t /D. The 
standard estimate is D t = ^u lms £, so D t /D = 
^Ref orc . For Run 3 we have Ref orc = 180, so we 
expect D t /D = 60. The actual value obtained 
by solving the passive scalar advection-diffusion 
equation simultaneously with Eqs. (l)-(3) is some- 
what smaller; see §3.4 where we find D t /D w 40 
for weak fields. This is probably due to the ab- 
sence of a proper inertial range. Ideally one would 
like to simulate higher levels of turbulence, which 
requires higher resolution. Certain questions can 
therefore not be addressed in a satisfactory man- 
ner, for example what are the spectral properties 
of the magnetic field, especially at large magnetic 
Prandtl numbers. Addressing this requires the 
presence of a sufficiently extended inertial range. 
Other aspects may very well be addressed, for ex- 
ample what is the behavior of the large scale field 
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Table 1 
Summary of runs. 





Run 1 


Run 21 


Run 2 


Run 3 


Run 3p 


Run 4 


Run 5 


Run 6 


Run 7 


mesh points 


30 3 


30 3 


go 3 


120 3 


120 3 


120 3 


120 3 


120 3 


120 3 




0.01 


0.005 


0.005 


0.002 


0.002 


0.02 


0.02 


0.001 


0.002 


V j 7) 


1 


1 


1 


1 


2 


20 


100 


1 


0.1 


k c 


5 






5 


5 


5 


5 


30 


5 




0.16 


0.23 


0.22 


0.29 




0.11 


0.114 


0.082 


0.29 


0.12 


0.15 


0.15 


0.18 


0.19 


0.10 


0.104 


0.062 


0.20 


Rm 


80 


200 


200 


600 


1200 


700 


3300 


400 


60 


-Rm,lin 


100 


300 


300 


900 


1800 


700 


3600 


500 


90 




3 


9 


9 


23 


46 


21 


112 


3 


2 


^?m,forc 


20 


60 


60 


180 


360 


140 


700 


17 


18 


^?m,forc,crit 


7.3 




6.9 


8.9 


8.9 


12 


12 




8.9 


A 


0.026 


0.06 


0.056 


0.067 




0.03 


0.04 


0.075 


0.036 


At 


0.24 


0.32 


0.34 


0.30 




0.34 


0.48 


0.19 


0.16 




0.18 


0.27 


0.28 


0.38 


0.40 


> 0.2 


> 0.21 


> 0.2 


0.18 




0.44 


0.75 


0.76 


1.27 


1.56 


0.7 


1.05 


1.5 


0.46 




0.80 


1.12 


1.12 


1.81 


1.81 


0.58 


0.58 


2.4 


1.76 


0.65 


0.78 


0.82 


1.23 


1.38 


0.55 


0.55 


1.8 


1.03 


(^ • «>lin 


0.13 


0.24 


0.24 


0.37 


0.37 


0.063 


0.063 


0.20 


0.36 


(w • u) sat 


0.08 


0.10 


0.11 


0.16 


0.17 


0.058 


0.055 


0.10 


0.19 


(J ■ B) max 


0.007 


0.025 


0.027 


0.07 


> 0.03 


0.025 


> 0.06 


0.040 


0.006 


(J ■ B> max /(J rms B rms ) 


0.22 


0.20 


0.20 


0.25 




0.23 


0.18 


0.25 


0.17 


|<J-B>/(B 2 )| lln 


0.9 


1.3 


1.1 


1-1.5 




2.2 


2.4 


5.5 


0.70 



questions it is important that there is sufficient 
scale separation between the energy carrying scale 
and the scale of the box. Furthermore it is im- 
portant to allow for sufficient separation between 
dynamic and resistive timescales in order to iden- 
tify properly the mechanisms affecting large scale 
dynamo action. A factor of 5 in scale separation 
seems to be a good compromise allowing still some 
degree of turbulent mixing to take place. 

3.1. The inverse cascade 

Consistent with previous studies in this field 
(e.g. Meneguzzi et al. 1981, Balsara & Pouquet 
1999), we find the development of large scale fields 
through an inverse cascade effect of the magnetic 
hclicity. This is best seen in the evolution of mag- 
netic energy spectra, £m(&); see Fig. 2. The ki- 
netic energy spectrum, EK(k), is also shown. 

The random initial condition has a k 4 power- 
spectrum, corresponding to a delta-correlated vec- 
tor potential. However, even though the initial 
field was smoothed, the spectrum is deformed sig- 
nificantly during the first few timesteps. During 
the interval 20 < t < 200 the spectrum is nearly 
shape invariant and grows at all scales at the same 
rate; see Fig. 2. This is typical of small-scale dy- 




1000 1200 



Fig. 1.— Evolution of magnetic energy and kinetic 
energy (per volume) for three values of the mag- 
netic Reynolds number, R m . Note that the growth 
rate and the saturation level increase with increas- 
ing R m . For reasons of clarity the curves of kinetic 
energy are not shown in the range 550 < t < 1000. 



and how does it depend on Reynolds and Prandtl 
numbers. We shall show that the spectral proper- 
ties are well converged at large scales, but the time 
scales for reaching a final state increase with mag- 
netic Reynolds number. In order to address these 
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Fig. 2. — Spectra of magnetic energy for Run 3 
during the initial growth phase at t = (dotted 
line), t = 20, 40, 220 (solid lines) and t = 
240 (dashed line). The kinetic energy spectrum 
(time averaged between 600 < t < 1000) is shown 
for comparison (dash-dotted line give the kinetic 
energy spectrum. 

namos (Kazantsev 1968). 

At t = 200 the magnetic energy approaches 
equipartition with the kinetic energy at small 
scales. After t = 240 the magnetic energy is 
in slight supercquipartition with the kinetic en- 
ergy at k > 10. This marks the beginning of a 
more complicated process (Fig. 3) during which 
the field at the largest possible scale (k = 1) 
continues to grow, but the field at intermediate 
wavenumbers (k — 2, 3, and 4) begins to decline. 
This process is essentially completed by the time 
t = 400. The significance of this process becomes 
clear when looking at the magnetic field evolution 
in real space. 

3.2. The emergence of a large scale field 

Although the magnetic field has reached equipar- 
tition already at t w 200, and its scale began to 
reach the largest possible scale of the box, it took 
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Fig. 3. — Spectra of magnetic energy for Run 3 
during the saturation phase at times between t = 
240 and 360. The time averaged spectrum for 
times between 600 and 1000 is shown as a thick 
line. (Only the range 1 < k < 10 is shown.) 

another 100 time units for the large scale field 
at scale k = 1 to fully develop and, more im- 
portantly, to suppress the power at intermediate 
scales. Looking at xy and yz cross-sections, two 
components of the field (B x and B z ) show the 
development of a large scale sinusoidal modula- 
tion through the entire box. In Fig. 4 we show xy 
slices of B x , but the yz cross-sections look qual- 
itatively similar, except for a 90° phase shift of 
B z in the y direction. This systematic phase shift 
is seen more clearly in a plot of the three field 
components averaged in the x and z directions; 
see Fig. 5. 

Although our forcing is isotropic, one particu- 
lar direction has been selected by the large scale 
magnetic field. In Runs 1 and 3 it was the y- 
direction, in Run 2 the z-direction, and in Run 5 
the x-direction. Which direction is selected de- 
pends on fine details of the initially random con- 
dition. Nevertheless, it is not until the time of 
saturation that the final selection is established, 
as can be seen in Fig. 6 where we plot the mag- 
netic energies of the mean field for the three possi- 
ble directions, denoted by E(K X ), E(K y ), E(K Z ). 
So, 



E(K X ) = 


((B) 2 yz ) x , 


(8) 


E(K y ) = 




(9) 


E(K Z ) = 


(( B )ly)z, 


(10) 



i=10 



(=100 





-3-2-10 1 2 3 
t=4-0Q 



-3 -2 -1 1 2 3 
(=1000 





Fig. 4. — Gray-scale images of cross-sections of 
B x (x, y,0) for Run 3 at different times showing 
the gradual build-up of the large scale magnetic 
field after t — 300. Dark (light) corresponds to 
negative (positive) values. Each image is scaled 
with respect to its min and max values. 



where the subscripts denote the directions of av- 
eraging. Thus, for any direction K, say the re- 
direction, we define corresponding mean fields by 
averaging in the two perpendicular directions (y 
and z in this case), and then we calculate their 
mean squared value. The time of selection, i.e. 
when one of the three E(Ki) becomes dominant, 
is earlier in the large Reynolds number cases. 

A quantity of theoretical interest is the ratio 
(B )/(B 2 ), which characterizes the fraction of 
space occupied by the large scale field. Initially 
this ratio is just ~ 2% (for Run 3) and ~ 0.7% (for 




Fig. 5. — The three magnetic field components av- 
eraged in the x and z directions at t = 1000. Note 
the 90° phase shift between {B x )(y) and (B z }(y), 
and that the functional form is nearly sinusoidal. 
Run 3. 

Run2 




450 



0.1000 
0.0100 
0.0010 



Run5 



10 



°" 4 



400 600 800 1000 1200 1400 1600 
t 



Fig. 6. — Magnetic energies (per volume) of those 
components of the large scale field whose wave vec- 
tors point in the x, y, or z direction. Which of the 
three possible force- free solutions is selected in the 
end depends on chance. 

Run 5), but later it begins to level off near 80% 
(Fig. 7). Most likely real astrophysical dynamos 
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are far less effective in producing such clean large 
scale fields, because in reality the helicity of the 
effective forcing will be far less than 100%. Nev- 
ertheless, it is important to notice that it is at 
least theoretically possible to achieve large scale 
field energies near or in excess of the kinetic en- 
ergy, even though the magnetic Reynolds number 
is reasonably high. 




200 400 600 800 1000 1200 
t 



Fig. 7.— Evolution of the ratio (B) 2 /(B 2 ), for 
Runs 3 and 5. Note that strong large scale fields 
are obtained even for large magnetic Reynolds and 
Prandtl numbers. 

We note that the phase of the large scale field 
may be drifting slowly as long as the large scale 
magnetic energy has not yet reached a fully steady 
state. In Run 3, for example, the phase was 
still drifting slowly in the y direction (speed ~ 
+1.5?7fci), but then it began to settle after t w 
1000. 

3.3. Spectral helicity and energy transfer 

The primary reason for the large scale field gen- 
eration is related to magnetic helicity conserva- 
tion. Once helicity is injected into the system, it 
tends to make the magnetic field also helical, as is 
seen from Fig. 8. For a closed or periodic system 
however the net magnetic helicity is conserved, ex- 
cept for diffusion at small scales, i.e. 

±(A.B)=-2 V (J.B). (11) 

Thus, if the magnetic field is to become helical, it 
must at first have equal amounts of positive and 
negative helicity. This feature, which is familiar 
in MHD (e.g. Seehafer 1996, Ji 1999), is also seen 
in hydrodynamical simulations (Biferale & Kerr 



1995). At later times, however, magnetic diffu- 
sion can destroy magnetic helicity at small scales, 
leaving magnetic helicity of opposite sign at large 
scales. This is best described by the evolution 
equation of the magnetic helicity spectrum which 
can be derived from the Fourier transformed in- 
duction equation (3), 

dAi 

-^ = e k - v j k , (12) 

where hats and subscripts k indicate three- 
dimensional Fourier transformation, and £ = 
u x B is the electromotive force. We write down 
the corresponding equation for the evolution of 
By, = ik x Ay. and derive from these equations 

the evolution equation for Ay, ■ By,, where aster- 
isks denote complex conjugation. Note that this 
is gauge invariant, because adding a gradient to 
Ay, corresponds to adding an ik ■ By, term which 
vanishes, because the magnetic field is solenoidal. 
We denote the real parts of the shell-integrated 
spectra of Ay, ■ By, and £y, ■ By, by HM(k,t) and 
Su(k,t), respectively, and obtain 

^H M (k,t) = 2S M (k,t)-2r]k 2 H M {k,t). (13) 

Note that J Hyi{k,t)dk = (A ■ B) and, because 
of helicity conservation, J SM(k,t)dk = 0, so it 
makes sense to write Eq. (13) in the form 

d d 

—H M (k,t) = -—G M (k,t) - 2rjk 2 H M (k,t), 

(14) 

where we have defined the spectral flux of helicity, 

G M (M) = / 2S M (k',t)dk', (15) 

which is plotted in Fig. 9 for different times. The 
magnetic helicity flux, Gyi(k,t), is found to be al- 
ways positive and its peak shifts from small scales 
(k rts 10) at early times to large scales (k = 2 — 3) 
at later times when the magnetic field becomes dy- 
namically important. Positive magnetic helicity is 
being produced on the right of the maximum of 
G-M.{k,t) and negative on the left; see Fig. 10. 
In view of the realizability condition, 

E M {k,t)>\kH M {k,t) (16) 

(e.g. Moffatt 1978), the spectral magnetic helic- 
ity can be viewed as the driver of spectral mag- 
netic energy: while small scale magnetic helicity 
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Fig. 8.— Evolution of —(A ■ B)/(B 2 }, compared 
with other magnetic length scales. Run 3. 
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Fig. 10. — Divergence of magnetic helicity flux, 
Sm (k), scaled by fc, so the area under the positive 
and negative parts of the curve are the same when 
plotted in a lin-log plot. Because of resistive losses 
of positive magnetic helicity at small scales, the 
resulting magnetic helicity, Hy[(k), is now dom- 
inated by negative helicity at large scales. The 
lower plot shows H M (k) , scaled by fc 3 to show the 
contributions from small scales. Both spectra are 
time averages over the interval 600 < t < 1000. 



Fig. 9. — Spectra of magnetic helicity flux at 
t = 10, 30, and 60 (dot-dashed), 100 an 130 (solid), 
160 and 200 (dashed), 300 (dotted), and for the 
time average between t = 600 and 1000 (thick 
line). 

is being destroyed, an equal amount gets into the 
large scales, and this must necessarily enhance the 
magnetic energy so as to satisfy (16). Indeed, in 
the present simulations the inequality (16) is al- 
most saturated at all scales, except at intermedi- 
ate scales 2 < k < 4; see Fig. 11. 

In order to determine the dominant interactions 
leading to the generation of the large scale field at 
k = 1 we now consider the spectral energy equa- 
tion, 

^-E M (k,t) = 2T M (k,t)-2r 1 k 2 E M (k,t), (17) 



where the transfer function of magnetic energy, 
T]yr(fc,i), is the shell-integrated spectrum of the 
real part of £^ ■ Jj*. Since £ = u x B, this corre- 
sponds really to a triple product, 

]T (u p xB q )-J* k , (18) 
k=p+q 

where the skew product can also be written as 
—UplJ^xBq), emphasizing that this term corre- 
sponds to the work done against the Lorentz force. 
In order to identify the dominant interactions we 
have calculated, in real space, the spectral transfer 
matrix 

T M (k,p,q,t) = -(u p -(J k xB q )}, (19) 

where angular brackets denote volume averages 
and the subscripts fc, p, and q denote Fourier fil- 
tering around the corresponding wavenumber (by 
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Fig. 11. — Time-averaged spectra of kinetic and 
magnetic energy, as well as kinetic and magnetic 
hclicity. Note that the magnetic energy exceeds 
the kinetic energy at k = 1, and that the inequal- 
ity (16) is almost saturated, except near k = 2 
and 3. The corresponding realizability condition 
for the kinetic helicity, on the other hand, is not 
very sharp. The dissipative cutoff wavenumber, 
(u; 2 /^ 2 ) 1 / 4 , is indicated by an arrow at the top. 
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Fig. 12. — Spectral energy transfer function 
T(k,p,t), normalized by (B 2 ) for three different 
times, for k = 1 and 2. Run 3. 

±1/2). (In this notation (Ak ■ Bfc), for example, 
is exactly the same as the helicity spectrum.) 
In Fig. 12 we show T(fc,p) = J2 q T(k,p, q), nor- 



malized by (.B 2 ) for the corresponding times, for 
k = 1 and 2. This function shows that most of the 
energy of the large scale field at k = 1 comes from 
velocity and magnetic field fluctuations at the forc- 
ing scale, k = 5. At early times this is also true 
of the energy of the magnetic field at k = 2, but 
at late times, t = 1000, the gain from the forcing 
scale, k = 5, has diminished, and instead there is 
now a net loss of energy into the next larger scale, 
k = 3, suggestive of a direct cascade operating at 
k = 2, and similarly at k = 3. 

The generation of large scale energy at k = 1 
through nonlocal inverse energy transfer is char- 
acteristic of the a-effect in mean-field electrody- 
namics. In the following we shall pursue this 
analogy further. It should be emphasized, how- 
ever, that without the simultaneous loss of en- 
ergy at the next smaller scales (here k — 2 and 
3) through direct energy transfer the k — 1 field 
would have been totally swamped by smaller scale 
fields. Thus, nonlinearity is quite crucial for this 
process to produce well defined large scale fields. 
Indeed, in the absence of the nonlinear term J x B 
in Eq. (2) the marked large scale pattern (Fig. 4) 
disappears within a turnover time. Recent numer- 
ical experiments have shown, however, that the 
ambipolar diffusion nonlinearity too leads to well 
defined large scale fields - even in the absence of 
the Lorentz force (Brandenburg & Subramanian 
2000). 

3.4. Mean-field interpretation 

In this subsection we adopt the hypothesis that 
the large scale component of the field at wavenum- 
ber k — 1 can be described in terms of mean-field 
theory. The magnetic field at other wavenumbcrs 
is (k > 2) is important for contributing to the a- 
effect and the turbulent magnetic diffusivity, r) t , 
but apart from that it is merely an extra source 
of noise as far as the dynamics of the large scale 
field is concerned. As we have seen in the previous 
section, this extra noise is automatically kept to a 
minimum due to direct cascade effects and transfer 
to kinetic energy during the saturation phase. 

According to mean-field theory for non-mirror 
symmetric isotropic homogeneous turbulence with 
no mean flow the mean magnetic field is governed 
by the equation 

^-B = aV xB + r lT V 2 B, (20) 



9 



where bars denote the mean fields, a and rjT = rj + 
% are constants, and r\t is the turbulent magnetic 
diffusivity. In general, these coefficients are not 
constant and depend for example on the magnetic 
field. (In our particular case the local magnetic 
energy density is however approximately uniform.) 
Furthermore, since the magnetic field is strong, a 
and T)t should really be replaced by tensors, but 
we shall ignore this additional modification except 
that we shall allow a and r) t to vary slowly in time 
as the magnetic field approaches saturation. This 
simplified form of nonlinearity may be justified by 
noting that the mean magnetic field looks nearly 
sinusoidal (Fig. 5). 

Equation (20) permits steady force-free solu- 
tions where the current helicity of the large scale 
field, J ■ B is given by (rfr/a)J . Apart from 
some common phase factor, the mean field de- 
picted in Fig. 5 is given by B(y) = (siny, 0, cosy), 
so J — (— sin y, 0, — cosy), corresponding to neg- 
ative helicity, and therefore a must be negative. 
This is in agreement with mean-field theory which 
predicts that a is a negative multiple of the resid- 
ual (kinetic minus current) helicity (e.g. Blackman 
& Chou 1997, Field, Blackman, & Chou 1999), 
which is positive in our case; see Table 1. 

If the wavevector of the large scale field is K y 
(as in the case discussed above), Eq. (20) becomes 

t x = +aB' z + Vt b' x , (21) 

B~ z = -aB~' x +r lT Bl, (22) 

where dots and primes denote differentiation with 
respect to t and y, respectively. Since a < 0, the 
solution can be written in the form 

B x {y,t) = b x {t)sm{y + 4>), (23) 

B z (y,t) =fe z (t)cos(y + 0), (24) 

where b x (t) and b z (t) are positive functions of time 
that satisfy 

b x = \a\b z - r] T b x , (25) 

b z = \a\b x -f] T b z . (26) 

In a steady state \a\ = rjT, and b x = b z . In order to 
find the actual values of \a\ and tjt during both the 
saturated steady state and the growth phase we 
can do a simple experiment: suppose we put b x = 



at some moment in time, then Eq. (25) would 
predict that b x starts to recover at the rate \a\b z , 
which allows us to estimate a. In practice we put 
b x = by subtracting B x from the x component 
of B at a certain time and restart the simulation 
with that field. 

In order to have a somewhat more precise esti- 
mate we need the solution to (25) and (26) for the 
initial condition b x (0) = 0; 

bx = e (\a\-vr)t _ g-dal+Tjr)^ ( 2 7) 

b x=e (M-rn)t + e -(\a\-hrr)t^ ( 2 8) 

where the amplitude is arbitrary in linear theory. 
Adding and subtracting (27) and (28) we can solve 
for | a | — t]t and \a\ + t]t, respectively. In terms of 
| a | and t]t separately, we have 

\ a \ = 2^ [+ H+b x + b z ) - H-b x + b z )} , (29) 

nr = ^[- H+bx + h) - H-b x + b z )] . (30) 

In practice we average the results of (29) and (30) 
over some 5-10 time units. We have applied this 
method to Runs 1-3 at times between t = 100 — 
600 and to Run 5 at times between t = 300 — 
1600 (Fig. 13). During these times the mean field 
was still evolving (Fig. 6), so at different times the 
mean magnetic field was different, which allows us 
to obtain the a(|B|) dependence. We take into 
account the fact that during the experiment the 
actual field is only l/\/2 of what it was before one 
of the two components of the mean field has been 
removed. We have then attempted a fit of the form 



l + a B (B )/B c 2 q 

where B% q = MoPo( M2 )sat- The result is shown in 
Fig. 14 and the coefficients ao and cub are listed in 
Table 2. One should note, however, that Eq. (31) 
does not accurately represent the actual data of 
Run 5. Nevertheless, it is clear that a-qucnching 
is enhanced for large values of R m , which may be 
described by a fit of the form 

"b w (i? m , forc /22) 1 - 35 . (32) 

Such a steep dependence of qb on i? mi forc was sug- 
gested by Vainshtein et al. (1993), although his 
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Fig. 14. — Results for a for different values of 
i? m using Eq. (29). The lines represent the fits 
described in the text. 



argument (see also Vainshtein & Cattaneo 1992) 
assumes the presence of strong small scale fluctu- 
ations (which is not the case here; see Fig. 7). 

In Fig. 15 we compare with the result for a ob- 
tained by just imposing a uniform magnetic field, 
B z, and calculating a simply as 



a" 



((« x B)) z /B . 



(33) 



Each point in Fig. 15 corresponds to a different 
run with given field strength B , but otherwise 
the same parameters as in Runs 1-3 and 5. This 
method was frequently used in the past (e.g. Bran- 
denburg et al. 1990, Tao et al. 1993, Cattaneo & 
Hughes 1996), but it is not a priori clear that one 
measures the same a as with the method explained 
above. Nevertheless, the two results appear to be 



qualitatively similar (cf. Figs 14 and 15) although 
there are some differences in the case where R m is 
very large (Table 2). The same values of «b are 
confirmed by yet another method that is explained 
below in S3. 6. 
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Fig. 15. — Results for a for different values of 
i? m using an imposed magnetic field. The lines 
represent the fits described in the text. 



Table 2: Coefficients for the a-quenching expres- 
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The only way a strongly i? m -depcndent a 
quenching can be compatible with the large scale 
field generation seen in the present simulations 
would be if rjt was also strongly quenched. [See 
Cattaneo & Vainshtein (1991) for two-dimensional 
simulations supporting the hypothesis of strong 7/ t 
quenching.] In Fig. 16 we compare the results 
obtained for 7/ T with the turbulent diffusion co- 
efficient for a passive scalar. The passive scalar 
diffusion coefficient is obtained by simultaneously 
solving an additional equation for the concentra- 
tion per unit mass, c, 



Dt p 



(34) 
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Fig. 16. — Results for t\t and Dt for different 
values of R m . The lines represent the fits described 
in the text. In the plot of ryr the asterisks denote 
|a| — A for the i? m! f rc = 120 run, which agrees 
reasonably well with t\t- 

where D = rj is chosen for Runs 1-3 and D = 2.5r/ 
for Run 5. The total (turbulent plus microscopic) 
passive scalar diffusion coefficient is obtained by 
measuring the rate at which a narrow gaussian 
distribution of c widens as time goes on. The re- 
sult is shown in the second panel of Fig. 16. Gen- 
erally the suppression of 77 t by magnetic fields is 
stronger than the suppression of D t . (For Run 3 
D T /D ranges between 40 for weak fields and 10 
for strong fields.) We note that a dependence of 
the form 

D T = D + D tQ /{l + D B (\B\)/B e(l ) (35) 

seems to fit at least the Dt data better than a 
quadratic dependence. Unlike the results for a, 
the dependence of Db on Reynolds number is here 
less strong. In the fits shown in Fig. 16 a good fit 
is -Db = 2 for all three runs, and Do = 0.02, 0.035, 
and 0.06 for Runs 1, 2, and 3, respectively. Run 5 
behaves differently, because in this large magnetic 
Prandtl number run there is no inertial range, and 
so we used D — 0.007 and D B = 6. 



Determining r/t from dynamo simulations is no- 
toriously difficult: rjt has to be determined si- 
multaneously with a from the electromotive force, 
where r/t multiplies a derivative of the field, so it 
is numerically more noisy. Nevertheless, one sees 
from Fig. 16 that % is quenched by more than a 
factor of ten. However, the functional form can- 
not be established from our data. Using for r]T a 
similar fit formula as Eq. (35) we have t]b = 16 
for all three runs, and 771-0 = 0.03, 0.07, and 0.2 
for Runs 1,2, and 3, respectively, and 77x0 =0.43, 
tjb = 60 for Run 5. However, for strong mag- 
netic fields, t]t levels off at a value of 0.01 (for 
Run 3) or 0.005 (for Run 5). These values are 
similar to the values of a. Thus, a ~ rixfci, which 
means that the mean-field dynamo has turned into 
a marginally critical state, which is indeed to be 
expected. As a consistency check for the directly 
obtained values of tjt we use the time dependent 
growth rate A(t) and show |a| — A (asterisks in 
Fig. 16) for Run 3. They agree reasonably well 
with t]t (except near (\B\) = 0.4£? cq ), suggest- 
ing that (25) and (26) are approximately satisfied 
with the coefficients obtained above. Later in §3.6 
we present a more accurate and self-consistent de- 
termination of the combined expression for a and 
?7t-quenching by fitting solutions of Eq. (20) to the 
actual evolution of the mean field. Those results 
support a quadratic quenching formula for both a 
and ?7t. 

3.5. Large scale separation 

In some earlier exploratory models we forced 
the flow at k = 10, which gave somewhat more 
room for the inverse cascade to develop, but 
less room for the direct cascade toward smaller 
wavenumbers. The latter means that the Reynolds 
number with respect to the forcing scale is smaller 
and the turbulent mixing properties, as quantified 
by the ratio D t /D, are poorer. Strong turbulent 
mixing was important in order to address the is- 
sue of Reynolds number dependent suppression of 
transport coefficients. In this subsection we shall 
accept poor mixing in favor of larger scale separa- 
tion and hence a significantly larger subrange in 
which the inverse cascade can develop. Thus, we 
now force the flow at k — 30 (Run 6). 

After saturation is reached, which happens first 
at some intermediate scale, there is a wave of spec- 
tral energy propagating toward smaller k. This is 
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similar to the results of closure models of Pouquet 
et al. (1976). However, unlike these models there 
does not seem to develop a fc -1 magnetic energy 
spectrum. Instead, there is only an envelope of 
the helicity wave that follows approximately a fc" 1 
power law (Fig. 17). As before, there is a strong 
built-up of magnetic energy at the largest scale, 
k = 1 , combined with a subsequent suppression of 
energy at intermediate scales (2 < k < 20). The 
ratio of the peak value at fci = 1 to the value of the 
secondary peak at the forcing scale kf = 30 scales 
like (fci/Zcf ) — 1 (Fig. 17). This is also the case in 
runs with kf — 5 (Fig. 11). 



i 




1 10 
k 



Fig. 17. — Magnetic energy spectra for Run 6 with 
forcing at k — 30. The times range from (dotted 
line) to 10, 30, 290 (solid lines). The thick solid 
line gives the final state at t = 1000. Note that at 
early times the spectra peaks at fc max ~ 7. The 
fc _1 and fc+ 3 / 2 slopes are given for orientation as 
dash-dotted lines. 

As a function of time, the spectral magnetic 
energy grows at the same rate at all values of k 
until saturation is reached. During early times the 
magnetic spectra peak at fc max s=y 7, which is also 
the wavenumber that reaches equipartition first, 
but then the field at this wavenumber decays to 
a somewhat smaller value while the contributions 



from the next smaller k grow and subsequently 
decay; see Fig. 18. Finally, when the energy at 
k = 1 reaches equipartition the energy of all larger 
ks becomes suppressed. 
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Fig. 18. — Evolution of spectral magnetic en- 
ergy for selected wavenumbers in a simulation with 
forcing at k — 30. 

Although a comparison with mean-field theory 
may be inappropriate, it is interesting to note that 
the existence of a fcmax is predicted from the theory 
of a 2 dynamos in infinite media (Moffatt 1978). 
The dispersion relation is 

X = \a\k - T] T k 2 , (36) 

where A is the growth rate. Maximum growth oc- 
curs at wavenumber 

fcmax = iMA?T, (37) 

and there the maximum growth rate is 

Amax = ^Ittlfcmax. (38) 

Since fc max and A max can be measured, we can de- 
termine 

\a\ = 2A max /|fc max , T] T = A max /|fc max (39) 

during the growth phase of the dynamo. With 
fcmax = 7 (Fig. 18) and A ma x = 0.07 (Fig. 17) we 
have |a| = 0.02 and rjr — 0.0014. In this case, the 
forcing occurs essentially in the dissipation range, 
so the turbulence has therefore poor mixing prop- 
erties {rj^/ri = 1.4). The value of a is about 7 
times smaller than that of Run 3 during the kine- 
matic regime, which seems reasonable. 
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In the nonlinear regime there are marked differ- 
ences between mean-field theory and simulations. 
Using a nonlinear mean-field model of the inverse 
cascade with a-effect, Galanti, Sulcm, & Gilbert 
(1991) found significant power at k > 2, which is 
in contrast to the behavior seen in Fig. 17. Con- 
sequently, the mean fields of Galanti, Gilbert, & 
Sulem (1990) look much less sinusoidal than in 
the simulations (Fig. 5). If one wanted to model 
this within the framework of an a-effect one would 
need to invoke a scale-dependent a integral kernel, 
where the dominant contributions to a come only 
from the largest scale of the system (Brandenburg 
& Sokoloff 2000). 

3.6. Reynolds number dependence and 
magnetic helicity 

In this section we first discuss to which ex- 
tent our results are affected by the limited reso- 
lution and finite magnetic Reynolds number. We 
then show that magnetic helicity conservation im- 
plies slow saturation of helical magnetic fields, and 
that this features is quantitatively reproduced by 
a mean-field model with i? m -dependent quenching 
of a and r\ t - 

In Fig. 19 we show energy spectra of the mag- 
netic field for three values of R m . The three curves 
differ by a certain factor, but are otherwise essen- 
tially the same at small wavenumbers. The main 
difference is, as expected, at the smaller scales: the 
spectra for large R m begin to show signs of an in- 
crtial range for values of k larger than the forcing 
wavenumber. 

Although the magnetic energy spectra of the 
statistically steady state seem to show conver- 
gence to a spectrum roughly compatible with 
fc~ 5 / 3 , there are some serious concerns about the 
timescale on which such a steady state is achieved. 
As was first pointed out by Berger (1984), in 
a closed box (periodic or perfectly conducting) 
there is an upper bound on the rate of change of 
the magnetic helicity. This is relevant, because the 
fields that are generated in the simulations have 
strong magnetic helicity: (A-B)/(B 2 ) w —0.8 /hi 
(Fig. 8). Most of the magnetic helicity is in the 
large scales (Fig. 11). In order to build up the he- 
licity at large scales we have to destroy magnetic 
helicity at small scales (Fig. 9). 

Open boundary conditions would help to get rid 
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Fig. 19. — Comparison of time averaged magnetic 
energy spectra for Runs 1-3 (t = 600 — 1000) 
and Run 5 (t = 1600). To compensate for dif- 
ferent field strengths and to make the spectra 
overlap at large scales, two of the three spectra 
have been multiplied by a scaling factor (x3.4 for 
Run 5, x2 for Run 2 and x5 for Run 1). There 
are signs of a gradual development of an inertial 
subrange for wavenumbers larger than the forc- 
ing scale. The fc -5 / 3 slope is shown for orienta- 
tion. The dissipative magnetic cutoff wavenum- 
bers, fed = (J 2 /?7 2 ) 1 / 4 , arc indicated by arrows at 
the top. 

of magnetic helicity (cf. Blackman & Field 2000), 
which may well be an important factor in more 
realistic simulations (e.g. Glatzmaier & Roberts 
1995, Brandenburg et al. 1995). However, the 
present simulations are for closed boxes and yet 
they do show strong fields, so we need to under- 
stand whether and how they have been affected by 
this constraint. 

The rate at which magnetic helicity can change, 
see Eq. (11), is constrained by the Schwarz in- 
equality, 



dt 



(A-B) 



2r)\(J-B}\ <2r/J rms B rms 



(40) 



(Berger 1984, Moffatt & Proctor 1985). In forced 
systems it is common (albeit not necessary) that 
the energy (-B 2 ms ) and dissipation rate (i] J r 2 ms ) arc 
independent of the Reynolds number. However, 
this would imply that the rate of change of (A ■ B) 
is limited resistively by a term proportional to 
77 1 / 2 . In our simulations, (J • B) / (J rms B ms ) has 
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a maximum of <~ 0.25 (see Table 1), so this limit 
does not seem to be saturated. Another estimate 
for the helicity dissipation is obtained by assuming 
that the magnetic spectrum has powerlaw behav- 
ior with Eyi{k) ~ k~ n for k < k^. The value 
of fc max follows by assuming that the Joule dis- 
sipation is independent of r\, which yields ~ 



i] 



-1/(3- 



■0 (~ ^-a/ 4 for n = 5/3). The helicity dis- 



,1/(3-7*) 



sipation is then proportional to <~ rj* 
which always tends to zero for small values of 
r\. For n — 5/3 the magnetic helicity dissipation 
scales like r/ 3 / 4 , which is faster than Berger's limit. 

There is a related and even more important 
consequence of helicity conservation. In a steady 
state, (A - B), which is gauge invariant and there- 
fore physically meaningful, must be constant, and 
hence (J-B) — ► 0. Thus, although we are going to 
generate a strong large scale field with significant 
current helicity in the large scales, the net cur- 
rent helicity must actually vanish. This can only 
happen if there is an equal amount of small scale 
current helicity of opposite sign, i.e. 

(j ■ b) w -(J ■ B) (in the steady state), (41) 

where angular brackets denote volume averages, 
bars denote the large scales at k = 1, and lower 
characters refer to contributions from all higher 
wavenumbers. The spectrum of the magnetic he- 
licity is k 2 times steeper than that of the cur- 
rent helicity, so (a • b) w (j ■ b)/kf, and therefore 

\(a-b)\«\(A-B)\- 

In the following we estimate the evolution of 

the energy density of the large scale field, B , us- 
ing the fact that in the present simulations a large 
fraction of the magnetic energy is contained in the 
large scale field (Fig. 7) and that magnetic helic- 
ity is strong (Fig. 8). This would be too ideal 
an assumption for astrophysical settings, but it is 
adequate for describing our present simulations. 
Thus, we may set 



k 1 (A-B)^k 1 (A-B) 



(J-B)/h. 

(42) 

For clarity we retain here the k\ factors, even 
though they are one. Using Eq. (11) this yields 



-2j]k\(B )+2 V k 1 \(j -b)\ 
which has the solution 

ft 

-2r)k(t / p lT\k\t' 



(43) 



(s 2 ) 



2rjkie 



j-b)\dt', (44) 



which is expected to apply after the time of sat- 
uration, t = i sa t, when (j ■ b) is approximately 
constant, so 



\(j-b)\ 



I _ e -2r;fc2(t- 



(45) 



To a good approximation we may assume (j ■ b) w 
Po(w ■ u). This means that the 'residual helicity' 
(Pouquet et al. 1976), (u> • u) — (j -b)/po, is small, 
which is indeed consistent with the present data 
and also with a results of Brandenburg & Subra- 
manian (2000), who used the ambipolar diffusion 
model of Subramanian (1999). The kinetic helic- 
ity can be approximated by (u> ■ u) s=s kf(u 2 ), and 
so the final field strength, B^ n , is given by 



(3-b)\ 



po- 



v Po{u) - = B eq -. 

(46) 

For Runs 1-3 and 5-7, the ratio of the actual val- 
ues of Bfi n to those anticipated from Eq. (46) are 
0.61, 0.76, 0.85, 1.16, 0.74, and 0.36. In the large 
R m cases (Runs 3 and 5) the ratio is close to one, 
confirming thus the assumption of small residual 
helicity in the saturated state. 

The resistively limited growth of B has imme- 
diate consequences for a. The evolution equation 
for (A ■ B) can be derived from Eq. (20), 

^-(A-B)=2a(B 2 )-2r lT (J-B). (47) 

In the strongly helical case considered here the 
magnetic helicity of the large scale field, (A-B), is 
very nearly equal to (A-B), so the right hand sides 
of (11) and (47) must be approximately equal, 
which leads to 



\a\-rrrh^ V \(J-B)\/(B >. 



(48) 



In order to check this relation we plot the two sides 
of (48) versus time (Fig. 20), where 7/t = V + Vt 
and | a | — r/ t ki has been obtained directly from 
£ = u x B assuming £ = aB — rjt J which, for 
fully helical mean fields, leads to 



\a\- Vt k 1 ^\(£-B)\/(B ). 



(49) 



(Note that here r/ t enters, not tjt- Recall also that 
in our case a < 0, so (£ ■ B) < 0.) The condi- 
tion (48) on a reflects the fact that the growth of 
the large scale field is limited by the microscopic 
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resistivity, as seen already from Eq. (45). Note, 
however, that this statement only applies to the 
present case of strong helicity and closed or peri- 
odic boxes. 
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Comparison of \a\ 



ifrki with rj(J • 

B)/(W) and rj cnp[-2r]kf(t - t sat )]. Note that 
| a | —r/^ki follows very closely the resistively dom- 
inated limit. Notice also the relatively large noise 
level in \a\ — t]t- Run 3. 

There have been previous attempts to incorpo- 
rate the conservation of helicity into the expression 
for a, which led to a dynamical dependence of a on 
time and field strength (Kleeorin, Rogachevskii, 
& Ruzmaikin 1995, Kleeorin et al. 2000). In the 
following we point out, however, that the helic- 
ity constraint, which leads to the prolonged sat- 
uration phase, is well described in terms of an 
a 2 -dynamo with resistively dominated quenching 
functions for both a and %, i.e. 

"o Vto 



a 



1 + a B B /B 2 



Vt = 



l + m B Z /B 2 



(50) 



where a B = f]B is assumed. Using the fact that 



—2 



the magnetic energy density of the mean field, B , 
is approximately uniform we can write Eq. (20) in 
the form 



A fie 2 
dt 



fig 2 ) = \ a o\^-' 
^ ' l + a B B 



ao\h - rjtokl —2 2 ^2 



+ a B B /B| q 
There is a steady solution with 



B -r]kfB . (51) 



where 



a B B 2 /B 2 cq = X/ V kl 



A = |a |fci - r] T0 kl 



(52) 
(53) 



is the kinematic growth rate of the dynamo and 
f7To = V + Vto- The time dependent equation (51) 
can be integrated to give the solution in the form 



B 2 /(l 



B 2 /B t 



2Xt 



(54) 



where £>;„; is the initial field strength and Bg n is 
the final field strength which is given by Eq. (52). 
The parameters Bi nl and Bg n can be obtained 
from Fig. 1. The kinematic growth rate A is the 
same for small and large scale fields and can hence 
be taken from Table 1. We emphasize that there 
is excellent agreement between the results of the 
simulation and the solution (54); see Fig. 21. 
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Fig. 21.— Evolution of (B 2 ) for Runs 1-3 and 5- 
7 (solid lines), compared with the solution (54) of 
the mean-field dynamo equations using resistively 
dominated a and r\t quenchings (dotted lines) 

Having obtained A and Ba n from the simulation 
we can use Eq. (52) to find for a B and rj B the 
result 

a B = t]b = -^2 (^) ■ ( 55 ) 



r]k{ 



Note that a B and r\ B are proportional to I/77, 
which implies i? m -dependent a and 7] t quenchings. 
The results for a B = r/ B are summarized in Ta- 
ble 3 for different runs. 
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To a good approximation B^ D /B 2 q w fcf/fci, see 
Eq. (46), so we expect as ~ X/(r]kik{). Since 
A w 0.3 x Wrms/^f we have 

0.3 U rms L /rn\ 
a B ~ T^o ~ °- 01 X #m,lin- (56) 

(27r)^ ?7 

Thus, as should scale with the magnetic Reynolds 
number based on the box scale - not the forc- 
ing scale. The difference is particularly evident 
when comparing with Run 6 where the forcing 
wavenumber is six times larger than in the other 
runs. For comparison we list in Table 3 also 
-Rm,iin/100 and i? m .f orc /20. The agreement be- 
tween Eqs. (55) and (56) is generally good, except 
for Run 7 where cub is larger than expected. This 
is probably because in this run with a small value 
of vjj] and small R m the value of Bfi n is smaller 
than expected from Eq. (46); see Eq. (55). 

In Table 3 we also give for completeness the 
values of B m i. Together with the values of r\ and 
A in Table 1 all the data entering Eq. (54) are now 
specified. It should be mentioned that the value 
of i?i n i is obtained from a fit and is only roughly 
comparable to the actual seed magnetic field in 
the simulation, where different initial structures 
are possible. 

Table 3: Coefficients for the a-quenching expres- 
sion obtained from Eq. (55) for runs with different 
values of magnetic Reynolds number. Note the 
rough agreement between the values of i? m j; n /100 
and «b for all runs except Run 7. 
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Run 3 


9 
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35 
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Run 6 
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11 


Run 7 
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2.8 


0.16 
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For weak fields Eq. (54) gives the usual ex- 
ponential growth, \B\ = Bi D i exp(At). For 
strong fields we recover Eq. (45) with t sat = 
A -1 ln(Sfi n /B in i) in the limit A > r\k\. 

3.7. Large magnetic Prandtl numbers 

The magnetic helicity constraint becomes more 
important as the magnetic Reynolds number is in- 
creased. So far we have mainly considered the case 



where u/r/ = 1. In the sun and many other as- 
trophysical bodies, vjr\ <C 1, but in the galaxy 
vjt] > 1 (e.g. Kulsrud & Andersen 1992). This 
may lead to a magnetic energy spectrum peak- 
ing at small scales (Kulsrud & Andersen 1992). 
However, although viscous damping will dissipate 
energy at small scales (Chandran 1998, Kinney 
et al. 1998), there is some recent evidence that 
the inverse cascade may no longer operate (Maron 
2000). These results have been obtained in the 
absence of net helicity. It will therefore be inter- 
esting to see whether in the presence of net helicity 
an inverse cascade is still possible when v/r]^> 1. 

In a preliminary attempt to clarify this ques- 
tion we have carried out simulations for vjr\ = 20 
(mainly by increasing the viscosity to v = 0.02; 
Run 4) and v/r] = 100 (where the magnetic dif- 
fusivity was lowered to r\ = 2 x 10~ 4 ; Run 5). 
The viscous cutoff wavenumber is then around 5, 
i.e. at the forcing scale, and the magnetic cutoff 
wavenumbers are 25 and 72, respectively. The res- 
olution for Run 5 may be insufficient and discreti- 
sation errors must play a role at small scales, but 
the images of the current density look reasonable 
(see below) and the evolution of the large scale 
field is also in agreement with expectations (§3.6). 

As seen from Fig. 22 the magnetic energy spec- 
trum does indeed peak at large wavenumbers ini- 
tially, although the magnetic energy spectrum 
does not scale like k 3 / 2 (Kulsrud & Andersen 
1992). Instead, the spectrum is close to fc 1 / 3 , 
which was also found during the kinematic stage 
of convective dynamos (Brandenburg et al. 1996). 
However, this result is inconclusive, because it 
could be an artifact of the lack of an inertial range 
in this run. In any case, at large wavenumbers 
the magnetic energy exceeds the kinetic energy, 
although at later times the kinetic energy is in- 
creased somewhat by magnetic forces. Especially 
at later times the magnetic energy is no longer 
dominated by small scales and the spectrum falls 
off more like fc~ 5 / 3 . The convergence to this pow- 
crlaw is evident when comparing Runs 3, 4, and 5 
(Fig. 23). Most importantly, there are now clear 
signs of an inverse cascade; see Fig. 24. 

In these runs with large magnetic Prandtl num- 
ber the current density shows strong filamentary 
structures that tend to be aligned with the local 
magnetic field direction, as seen in Fig. 25. The 
resulting anisotropy affects particularly the small 
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Fig. 22. — Spectra of magnetic and kinetic energy 
and helicity for a large magnetic Prandtl number 
at early and late times, v = 0.02, f) = 2x 10~ 4 . 
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Fig. 23. — Comparison of the magnetic energy 
spectra for Runs 3, 4 and 5. The magnetic cutoff 
wavenumbers are 25, 26, and 72, as indicated by 
arrows at the top. 

scales (Goldreich & Sridhar 1997, Maron 2000). 
Note that this type of anisotropy cannot be cap- 
tured by closure models (e.g. Pouquet et al. 1976). 



Fig. 24. — Images of the B y and B z components 
of the magnetic field in an arbitrarily chosen xy 
plane. Run 5, t = 1600. 



2.0 x 0,0 x 2.0 




Fig. 25. — Images of the current density, J Xl in the 
plane of the large scale magnetic field. The field 
direction is shown as a vector. Run 5, t = 1600. 

The main shortcoming of the present large 
Prandtl number calculations is that the viscous 
dissipation cutoff wavenumber is so small that it 
lies in the range of the forcing scale, so no iner- 
tial range in the kinetic energy is possible. At the 
same time, of course, the range of scales available 
to the magnetic field is still not large enough to 
establish a fc 3 / 2 scaling at early times. 

4. Conclusions 

The main conclusion to be drawn from this 
work is that in the presence of net magnetic helic- 
ity there is a gradual built-up of a nearly force- 
free magnetic field at the largest possible scale 
of the system. In our periodic calculations this 
corresponds to a sinusoidal one-dimensional Bel- 
trami field, e.g. B ~ (cos z, sin z, 0), which is of 
course locally strongly distorted by the turbulence. 
Nevertheless, the presence of the large scale field 
is clearly seen without averaging (Fig. 4). We 
emphasize that this result is numerically robust: 
the relative dominance of magnetic energy at the 
smallest wavenumber is independent of resolution 
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(Fig. 19) and independent of the degree of scale 
separation (Fig. 17). So, the effect is seen equally 
well at resolutions ranging from 30 3 to 120 3 mesh- 
points, and at forcing wavenumbers ranging from 5 
to 30. However, the time it takes to establish such 
large scale fields increases with the ohmic diffusion 
time. We also note that the results are not very 
sensitive to the choice of the forcing function: a 
forcing function that is nearly delta-correlated in 
space, but still strongly helical, yields very similar 
results. In the absence of net helicity, however, 
no large scale field is generated. Likewise, if the 
forcing is made non-helical the large scale field dis- 
appears. 

An important property of the turbulence is that 
once the large scale field is established, it can sup- 
press magnetic energy on scales smaller than the 
largest one. This leads to something like a 'self- 
cleaning' processes. This is also seen in histograms 
of the magnetic field which are, for the present 
simulations, more nearly gaussian (with one hump 
perpendicular to the field, and two in the direction 
of the field). This is similar to other simulations 
with large scale dynamo action (see Brandenburg 
et al. 1995), but very different from simulations of 
small scale dynamo action where the histograms of 
the field components show stretched exponentials 
(Brandenburg et al. 1996), which can also be seen 
in the present simulations, but only at early times. 

Our simulations show that most of the energy 
input to the large scale field comes from small 
scales. This type of nonlocal spectral energy trans- 
fer is suggestive of an a-effect that could be re- 
sponsible for the field generation, rather than a lo- 
cal inverse cascade, which transports energy from 
k = 2 to k = 1, for example. Although a local 
inverse cascade seems to occur at early times, i.e. 
before the magnetic field is fully established, once 
the field is strong the magnetic energy at k — 2 is 
actually cascaded to k — 3 and/or transferred to 
kinetic energy, both of which are probably impor- 
tant for the 'self-cleaning' processes. A sketch of 
the anticipated energy transfer properties is given 
in Fig. 26. 

We point out that the present simulations must 
not be regarded as local in the sense of represent- 
ing only a small chunk of a larger system, be- 
cause the field structure depends crucially on the 
size of the box. Instead, they should be viewed 
as global within the geometry considered. With 
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Fig. 26. — Sketch illustrating direct and inverse 
cascade processes in helical MHD turbulence. 

other boundary conditions or in different geome- 
tries the shape of the large scale field will be dif- 
ferent. In the case of a sphere, for example, no 
perfectly force-free field is possible, but the field 
may be nearly force-free. An example may be the 
field obtained in hydromagnetic calculations with 
a-effect (Proctor 1977), where field saturation oc- 
curs through the Lorentz force of the large scale 
field. In these calculations the magnetic saturation 
field strength is relatively large, which reflects the 
fact that the field is indeed nearly force- free. 

It should be emphasized that the overall growth 
of the large scale field and the saturation phase 
of the dynamo are well described by a simple 
a 2 -dynamo with a and % coefficients that arc 
quenched in a i? m -dependent fashion; see Eqs. (50) 
and (55). The reason such a dynamo can still sat- 
urate is because of the presence of microscopic dif- 
fusion, and it is this what causes the saturation to 
happen so slowly. The excellent agreement in the 
evolution toward saturation between both the sim- 
ulation and the mean-field model is an indication 
that the simple quadratic quenching formula is ac- 
tually correct. For example a cubic nonlincarity 
(Moffatt 1972, Riidiger 1974) would lead to differ- 
ent behavior and would not have the correct resis- 
tive relaxation asymptotics consistent with helic- 
ity conservation (Brandenburg 2000). 

The slow resistive field evolution past equiparti- 
tion has become particularly clear in Run 5, where 
the final selection of the large scale field structure 
occurred rather late (after t w 1200, corresponding 
to about 100 turnover times; Fig. 6). By contrast, 
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in Run 3, where the magnetic Reynolds number 
was about 6 times smaller, the large scale field was 
fully developed by the time t w 400, correspond- 
ing to about 50 turnover times. In stars the typical 
magnetic Reynolds numbers are at least another 
six orders of magnitude larger than in Run 5, so a 
large scale field, if generated by an a-effect, would 
require ~ 10 8 turnover times or ~ 3 x 10 6 yr (as- 
suming a turnover time of 10 days). In the case of 
the sun this estimate would be reduced by another 
factor of 100 (Brandenburg et al. 2000), because 
differential rotation contributes to non-helical field 
generation, so the resulting fields are only par- 
tially subject to the helicity constraint. Since even 
the youngest protostars are older than 3 x 10 4 yr 
the aw-dynamo may still be responsible for field 
generation in these bodies. For galaxies, on the 
other hand, the magnetic Reynolds numbers are 
by another seven orders of magnitude larger than 
in stars, making here the case for an aw-dynamo 
more doubtful. However, this assumes that the 
conclusions from models with closed or periodic 
boundaries apply to galaxies, and that the micro- 
scopic resistivity is not enhanced during reconnec- 
tion [see Ji et al. (1998) for anomalous resistivities 
in a laboratory reconncction experiment] . 

There is now also some evidence that in os- 
cillatory dynamos of aw-type the cycle period is 
not strongly affected by the helicity timescale con- 
straint (Brandenburg et al. 2000). This could 
be related to the fact that with shear the large 
scale field is no longer fully force-free and that in 
that case the turbulent magnetic diffusivity is only 
partially suppressed (Gruzinov & Diamond 1996). 
However, the case for aw-dynamo action in stars, 
galaxies or accretion discs is by no means settled. 
Firstly, proposals have been made for nonhelical 
large scale dynamo action (Vishniac & Cho 2000, 
Zheligovsky, Podvigina, & Frisch 2000), which 
may avoid the problems that aw-dynamos have. 
Secondly, real astrophysical bodies do have open 
boundaries and may get rid of small scale helicity 
rather rapidly (Berger & Ruzmaikin 2000). Indi- 
cations are, however, that open boundaries also 
produce significant losses at large scales, which 
lowers the overall dynamo efficiency; preliminary 
results are reported in Brandenburg (2000). 
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